library_load <- suppressMessages(
list(
# Seurat
library(Seurat),
# Data
library(tidyverse),
# Reticulate
library(reticulate),
# Plotting
library(ggplot2)
)
)
random_seed <- 42
set.seed(random_seed)
options(warn=-1)
# Set working directory to project root
setwd("/research/peer/fdeckert/FD20200109SPLENO")
# Source files
source("plotting_global.R")
source("bin/SeuratQC.R")
# SingleRPlot
source("bin/SingleRQC.R")
so <- readRDS("data/object/qc.rds")
so$integrate <- paste0(so$treatment, "-", so$sample_rep)
so <- SplitObject(so, split.by="integrate")
so <- so[c("NaCl-Rep2", "NaCl-Rep1", "CpG-Rep2", "CpG-Rep1")]
so <- lapply(so, function(so) {
so <- SCTransform(so, return.only.var.genes=FALSE, variable.features.n=8000, verbose=FALSE)
so <- RunPCA(so, npcs=100, assay="SCT", verbose=FALSE)
so <- FindNeighbors(so, dims=1:30, k.param=50, verbose=FALSE)
so <- FindClusters(so, verbose=FALSE, resolution=1, algorithm=1, group.singletons=TRUE)
so <- RunUMAP(so, dims=1:100, n.neighbors=50, verbose=FALSE)
return(so)
}
)
saveRDS(so, "data/object/so_sct.rds")
# so <- readRDS("data/object/so_sct.rds")
options(repr.plot.width=20, repr.plot.height=15)
dplot_1 <- dplot(so[[1]], reduction="umap", group_by="seurat_clusters") + ggtitle("Louvain")
dplot_2 <- dplot(so[[1]], reduction="umap", group_by="treatment") + scale_color_manual(values=color$treatment) + ggtitle("Treatment")
dplot_3 <- dplot(so[[1]], reduction="umap", group_by="sample_rep") + ggtitle("Replicate")
dplot_4 <- dplot(so[[1]], reduction="umap", group_by="tissue") + scale_color_manual(values=color$tissue) + ggtitle("Tissue")
dplot_5 <- dplot(so[[1]], reduction="umap", group_by="cc_phase_class") + scale_color_manual(values=color$cc_phase_class) + ggtitle("Cell cycle")
dplot_6 <- dplot(so[[1]], reduction="umap", group_by="label_fine_haemosphere") + scale_color_manual(values=color$label_fine_haemosphere) + ggtitle("Haemosphere label")
fplot_1 <- fplot(so[[1]], reduction="umap", features="pMt_RNA") + scale_colour_viridis_c(option="inferno") + ggtitle("Percentage Mt")
fplot_2 <- fplot(so[[1]], reduction="umap", features="pHb_RNA") + scale_colour_viridis_c(option="inferno") + ggtitle("Percentage Hb")
fplot_3 <- fplot(so[[1]], reduction="umap", features="pRb_RNA") + scale_colour_viridis_c(option="inferno") + ggtitle("Percentage Rb")
fplot_4 <- fplot(so[[1]], reduction="umap", features="nCount_RNA") + scale_colour_viridis_c(option="inferno") + ggtitle("UMI")
fplot_5 <- fplot(so[[1]], reduction="umap", features="nFeature_RNA") + scale_colour_viridis_c(option="inferno") + ggtitle("Feature")
dplot_comb <- dplot_1 + dplot_2 + dplot_3 + dplot_4 + dplot_5 + dplot_6 + fplot_1 + fplot_2 + fplot_3 + fplot_4 + fplot_5 + plot_layout(ncol=4)
dplot_comb
Scale for 'colour' is already present. Adding another scale for 'colour', which will replace the existing scale. Scale for 'colour' is already present. Adding another scale for 'colour', which will replace the existing scale. Scale for 'colour' is already present. Adding another scale for 'colour', which will replace the existing scale. Scale for 'colour' is already present. Adding another scale for 'colour', which will replace the existing scale. Scale for 'colour' is already present. Adding another scale for 'colour', which will replace the existing scale.
options(repr.plot.width=20, repr.plot.height=15)
dplot_1 <- dplot(so[[2]], reduction="umap", group_by="seurat_clusters") + ggtitle("Louvain")
dplot_2 <- dplot(so[[2]], reduction="umap", group_by="treatment") + scale_color_manual(values=color$treatment) + ggtitle("Treatment")
dplot_3 <- dplot(so[[2]], reduction="umap", group_by="sample_rep") + ggtitle("Replicate")
dplot_4 <- dplot(so[[2]], reduction="umap", group_by="tissue") + scale_color_manual(values=color$tissue) + ggtitle("Tissue")
dplot_5 <- dplot(so[[2]], reduction="umap", group_by="cc_phase_class") + scale_color_manual(values=color$cc_phase_class) + ggtitle("Cell cycle")
dplot_6 <- dplot(so[[2]], reduction="umap", group_by="label_fine_haemosphere") + scale_color_manual(values=color$label_fine_haemosphere) + ggtitle("Haemosphere label")
fplot_1 <- fplot(so[[2]], reduction="umap", features="pMt_RNA") + scale_colour_viridis_c(option="inferno") + ggtitle("Percentage Mt")
fplot_2 <- fplot(so[[2]], reduction="umap", features="pHb_RNA") + scale_colour_viridis_c(option="inferno") + ggtitle("Percentage Hb")
fplot_3 <- fplot(so[[2]], reduction="umap", features="pRb_RNA") + scale_colour_viridis_c(option="inferno") + ggtitle("Percentage Rb")
fplot_4 <- fplot(so[[2]], reduction="umap", features="nCount_RNA") + scale_colour_viridis_c(option="inferno") + ggtitle("UMI")
fplot_5 <- fplot(so[[2]], reduction="umap", features="nFeature_RNA") + scale_colour_viridis_c(option="inferno") + ggtitle("Feature")
dplot_comb <- dplot_1 + dplot_2 + dplot_3 + dplot_4 + dplot_5 + dplot_6 + fplot_1 + fplot_2 + fplot_3 + fplot_4 + fplot_5 + plot_layout(ncol=4)
dplot_comb
Scale for 'colour' is already present. Adding another scale for 'colour', which will replace the existing scale. Scale for 'colour' is already present. Adding another scale for 'colour', which will replace the existing scale. Scale for 'colour' is already present. Adding another scale for 'colour', which will replace the existing scale. Scale for 'colour' is already present. Adding another scale for 'colour', which will replace the existing scale. Scale for 'colour' is already present. Adding another scale for 'colour', which will replace the existing scale.
options(repr.plot.width=20, repr.plot.height=15)
dplot_1 <- dplot(so[[3]], reduction="umap", group_by="seurat_clusters") + ggtitle("Louvain")
dplot_2 <- dplot(so[[3]], reduction="umap", group_by="treatment") + scale_color_manual(values=color$treatment) + ggtitle("Treatment")
dplot_3 <- dplot(so[[3]], reduction="umap", group_by="sample_rep") + ggtitle("Replicate")
dplot_4 <- dplot(so[[3]], reduction="umap", group_by="tissue") + scale_color_manual(values=color$tissue) + ggtitle("Tissue")
dplot_5 <- dplot(so[[3]], reduction="umap", group_by="cc_phase_class") + scale_color_manual(values=color$cc_phase_class) + ggtitle("Cell cycle")
dplot_6 <- dplot(so[[3]], reduction="umap", group_by="label_fine_haemosphere") + scale_color_manual(values=color$label_fine_haemosphere) + ggtitle("Haemosphere label")
fplot_1 <- fplot(so[[3]], reduction="umap", features="pMt_RNA") + scale_colour_viridis_c(option="inferno") + ggtitle("Percentage Mt")
fplot_2 <- fplot(so[[3]], reduction="umap", features="pHb_RNA") + scale_colour_viridis_c(option="inferno") + ggtitle("Percentage Hb")
fplot_3 <- fplot(so[[3]], reduction="umap", features="pRb_RNA") + scale_colour_viridis_c(option="inferno") + ggtitle("Percentage Rb")
fplot_4 <- fplot(so[[3]], reduction="umap", features="nCount_RNA") + scale_colour_viridis_c(option="inferno") + ggtitle("UMI")
fplot_5 <- fplot(so[[3]], reduction="umap", features="nFeature_RNA") + scale_colour_viridis_c(option="inferno") + ggtitle("Feature")
dplot_comb <- dplot_1 + dplot_2 + dplot_3 + dplot_4 + dplot_5 + dplot_6 + fplot_1 + fplot_2 + fplot_3 + fplot_4 + fplot_5 + plot_layout(ncol=4)
dplot_comb
Scale for 'colour' is already present. Adding another scale for 'colour', which will replace the existing scale. Scale for 'colour' is already present. Adding another scale for 'colour', which will replace the existing scale. Scale for 'colour' is already present. Adding another scale for 'colour', which will replace the existing scale. Scale for 'colour' is already present. Adding another scale for 'colour', which will replace the existing scale. Scale for 'colour' is already present. Adding another scale for 'colour', which will replace the existing scale.
options(repr.plot.width=20, repr.plot.height=15)
dplot_1 <- dplot(so[[4]], reduction="umap", group_by="seurat_clusters") + ggtitle("Louvain")
dplot_2 <- dplot(so[[4]], reduction="umap", group_by="treatment") + scale_color_manual(values=color$treatment) + ggtitle("Treatment")
dplot_3 <- dplot(so[[4]], reduction="umap", group_by="sample_rep") + ggtitle("Replicate")
dplot_4 <- dplot(so[[4]], reduction="umap", group_by="tissue") + scale_color_manual(values=color$tissue) + ggtitle("Tissue")
dplot_5 <- dplot(so[[4]], reduction="umap", group_by="cc_phase_class") + scale_color_manual(values=color$cc_phase_class) + ggtitle("Cell cycle")
dplot_6 <- dplot(so[[4]], reduction="umap", group_by="label_fine_haemosphere") + scale_color_manual(values=color$label_fine_haemosphere) + ggtitle("Haemosphere label")
fplot_1 <- fplot(so[[4]], reduction="umap", features="pMt_RNA") + scale_colour_viridis_c(option="inferno") + ggtitle("Percentage Mt")
fplot_2 <- fplot(so[[4]], reduction="umap", features="pHb_RNA") + scale_colour_viridis_c(option="inferno") + ggtitle("Percentage Hb")
fplot_3 <- fplot(so[[4]], reduction="umap", features="pRb_RNA") + scale_colour_viridis_c(option="inferno") + ggtitle("Percentage Rb")
fplot_4 <- fplot(so[[4]], reduction="umap", features="nCount_RNA") + scale_colour_viridis_c(option="inferno") + ggtitle("UMI")
fplot_5 <- fplot(so[[4]], reduction="umap", features="nFeature_RNA") + scale_colour_viridis_c(option="inferno") + ggtitle("Feature")
dplot_comb <- dplot_1 + dplot_2 + dplot_3 + dplot_4 + dplot_5 + dplot_6 + fplot_1 + fplot_2 + fplot_3 + fplot_4 + fplot_5 + plot_layout(ncol=4)
dplot_comb
Scale for 'colour' is already present. Adding another scale for 'colour', which will replace the existing scale. Scale for 'colour' is already present. Adding another scale for 'colour', which will replace the existing scale. Scale for 'colour' is already present. Adding another scale for 'colour', which will replace the existing scale. Scale for 'colour' is already present. Adding another scale for 'colour', which will replace the existing scale. Scale for 'colour' is already present. Adding another scale for 'colour', which will replace the existing scale.
features_2000 <- SelectIntegrationFeatures(object.list=so, nfeatures=2000)
so <- PrepSCTIntegration(object.list=so, anchor.features=features_2000)
anchors <- FindIntegrationAnchors(object.list=so, normalization.method="SCT", anchor.features=features_2000, verbose=FALSE)
so_int <- IntegrateData(anchorset=anchors, normalization.method="SCT", verbose=FALSE)
so_int <- RunPCA(so_int, npcs=100, assay="integrated", verbose=FALSE)
so_int <- FindNeighbors(so_int, dims=1:30, k.param=50, verbose=FALSE)
so_int <- FindClusters(so_int, verbose=FALSE, resolution=1, algorithm=1, group.singletons=TRUE)
so_int <- RunUMAP(so_int, dims=1:100, n.neighbors=50, verbose=FALSE)
saveRDS(so_int, "data/object/so_sct_int_hvg2000.rds")
# so_int <- readRDS("data/object/so_sct_int_hvg2000.rds")
options(repr.plot.width=20, repr.plot.height=15)
dplot_1 <- dplot(so_int, reduction="umap", group_by="seurat_clusters") + ggtitle("Louvain")
dplot_2 <- dplot(so_int, reduction="umap", group_by="treatment") + scale_color_manual(values=color$treatment) + ggtitle("Treatment")
dplot_3 <- dplot(so_int, reduction="umap", group_by="sample_rep") + ggtitle("Replicate")
dplot_4 <- dplot(so_int, reduction="umap", group_by="tissue") + scale_color_manual(values=color$tissue) + ggtitle("Tissue")
dplot_5 <- dplot(so_int, reduction="umap", group_by="cc_phase_class") + scale_color_manual(values=color$cc_phase_class) + ggtitle("Cell cycle")
dplot_6 <- dplot(so_int, reduction="umap", group_by="label_fine_haemosphere") + scale_color_manual(values=color$label_fine_haemosphere) + ggtitle("Haemosphere label")
fplot_1 <- fplot(so_int, reduction="umap", features="pMt_RNA") + scale_colour_viridis_c(option="inferno") + ggtitle("Percentage Mt")
fplot_2 <- fplot(so_int, reduction="umap", features="pHb_RNA") + scale_colour_viridis_c(option="inferno") + ggtitle("Percentage Hb")
fplot_3 <- fplot(so_int, reduction="umap", features="pRb_RNA") + scale_colour_viridis_c(option="inferno") + ggtitle("Percentage Rb")
fplot_4 <- fplot(so_int, reduction="umap", features="nCount_RNA") + scale_colour_viridis_c(option="inferno") + ggtitle("UMI")
fplot_5 <- fplot(so_int, reduction="umap", features="nFeature_RNA") + scale_colour_viridis_c(option="inferno") + ggtitle("Feature")
dplot_comb <- dplot_1 + dplot_2 + dplot_3 + dplot_4 + dplot_5 + dplot_6 + fplot_1 + fplot_2 + fplot_3 + fplot_4 + fplot_5 + plot_layout(ncol=4)
dplot_comb
Scale for 'colour' is already present. Adding another scale for 'colour', which will replace the existing scale. Scale for 'colour' is already present. Adding another scale for 'colour', which will replace the existing scale. Scale for 'colour' is already present. Adding another scale for 'colour', which will replace the existing scale. Scale for 'colour' is already present. Adding another scale for 'colour', which will replace the existing scale. Scale for 'colour' is already present. Adding another scale for 'colour', which will replace the existing scale.
features_4000 <- SelectIntegrationFeatures(object.list=so, nfeatures=4000)
so <- PrepSCTIntegration(object.list=so, anchor.features=features_4000)
anchors <- FindIntegrationAnchors(object.list=so, normalization.method="SCT", anchor.features=features_4000, verbose=FALSE)
so_int <- IntegrateData(anchorset=anchors, normalization.method="SCT", verbose=FALSE)
so_int <- RunPCA(so_int, npcs=100, assay="integrated", verbose=FALSE)
so_int <- FindNeighbors(so_int, dims=1:30, k.param=50, verbose=FALSE)
so_int <- FindClusters(so_int, verbose=FALSE, resolution=1, algorithm=1, group.singletons=TRUE)
so_int <- RunUMAP(so_int, dims=1:100, n.neighbors=50, verbose=FALSE)
saveRDS(so_int, "data/object/so_sct_int_hvg4000.rds")
# so_int <- readRDS("data/object/so_sct_int_hvg4000.rds")
options(repr.plot.width=20, repr.plot.height=15)
dplot_1 <- dplot(so_int, reduction="umap", group_by="seurat_clusters") + ggtitle("Louvain")
dplot_2 <- dplot(so_int, reduction="umap", group_by="treatment") + scale_color_manual(values=color$treatment) + ggtitle("Treatment")
dplot_3 <- dplot(so_int, reduction="umap", group_by="sample_rep") + ggtitle("Replicate")
dplot_4 <- dplot(so_int, reduction="umap", group_by="tissue") + scale_color_manual(values=color$tissue) + ggtitle("Tissue")
dplot_5 <- dplot(so_int, reduction="umap", group_by="cc_phase_class") + scale_color_manual(values=color$cc_phase_class) + ggtitle("Cell cycle")
dplot_6 <- dplot(so_int, reduction="umap", group_by="label_fine_haemosphere") + scale_color_manual(values=color$label_fine_haemosphere) + ggtitle("Haemosphere label")
fplot_1 <- fplot(so_int, reduction="umap", features="pMt_RNA") + scale_colour_viridis_c(option="inferno") + ggtitle("Percentage Mt")
fplot_2 <- fplot(so_int, reduction="umap", features="pHb_RNA") + scale_colour_viridis_c(option="inferno") + ggtitle("Percentage Hb")
fplot_3 <- fplot(so_int, reduction="umap", features="pRb_RNA") + scale_colour_viridis_c(option="inferno") + ggtitle("Percentage Rb")
fplot_4 <- fplot(so_int, reduction="umap", features="nCount_RNA") + scale_colour_viridis_c(option="inferno") + ggtitle("UMI")
fplot_5 <- fplot(so_int, reduction="umap", features="nFeature_RNA") + scale_colour_viridis_c(option="inferno") + ggtitle("Feature")
dplot_comb <- dplot_1 + dplot_2 + dplot_3 + dplot_4 + dplot_5 + dplot_6 + fplot_1 + fplot_2 + fplot_3 + fplot_4 + fplot_5 + plot_layout(ncol=4)
dplot_comb
Scale for 'colour' is already present. Adding another scale for 'colour', which will replace the existing scale. Scale for 'colour' is already present. Adding another scale for 'colour', which will replace the existing scale. Scale for 'colour' is already present. Adding another scale for 'colour', which will replace the existing scale. Scale for 'colour' is already present. Adding another scale for 'colour', which will replace the existing scale. Scale for 'colour' is already present. Adding another scale for 'colour', which will replace the existing scale.
features_6000 <- SelectIntegrationFeatures(object.list=so, nfeatures=6000)
so <- PrepSCTIntegration(object.list=so, anchor.features=features_6000)
anchors <- FindIntegrationAnchors(object.list=so, normalization.method="SCT", anchor.features=features_6000, verbose=FALSE)
so_int <- IntegrateData(anchorset=anchors, normalization.method="SCT", verbose=FALSE)
so_int <- RunPCA(so_int, npcs=100, assay="integrated", verbose=FALSE)
so_int <- FindNeighbors(so_int, dims=1:30, k.param=50, verbose=FALSE)
so_int <- FindClusters(so_int, verbose=FALSE, resolution=1, algorithm=1, group.singletons=TRUE)
so_int <- RunUMAP(so_int, dims=1:100, n.neighbors=50, verbose=FALSE)
saveRDS(so_int, "data/object/so_sct_int_hvg6000.rds")
# so_int <- readRDS("data/object/so_sct_int_hvg6000.rds")
options(repr.plot.width=20, repr.plot.height=15)
dplot_1 <- dplot(so_int, reduction="umap", group_by="seurat_clusters") + ggtitle("Louvain")
dplot_2 <- dplot(so_int, reduction="umap", group_by="treatment") + scale_color_manual(values=color$treatment) + ggtitle("Treatment")
dplot_3 <- dplot(so_int, reduction="umap", group_by="sample_rep") + ggtitle("Replicate")
dplot_4 <- dplot(so_int, reduction="umap", group_by="tissue") + scale_color_manual(values=color$tissue) + ggtitle("Tissue")
dplot_5 <- dplot(so_int, reduction="umap", group_by="cc_phase_class") + scale_color_manual(values=color$cc_phase_class) + ggtitle("Cell cycle")
dplot_6 <- dplot(so_int, reduction="umap", group_by="label_fine_haemosphere") + scale_color_manual(values=color$label_fine_haemosphere) + ggtitle("Haemosphere label")
fplot_1 <- fplot(so_int, reduction="umap", features="pMt_RNA") + scale_colour_viridis_c(option="inferno") + ggtitle("Percentage Mt")
fplot_2 <- fplot(so_int, reduction="umap", features="pHb_RNA") + scale_colour_viridis_c(option="inferno") + ggtitle("Percentage Hb")
fplot_3 <- fplot(so_int, reduction="umap", features="pRb_RNA") + scale_colour_viridis_c(option="inferno") + ggtitle("Percentage Rb")
fplot_4 <- fplot(so_int, reduction="umap", features="nCount_RNA") + scale_colour_viridis_c(option="inferno") + ggtitle("UMI")
fplot_5 <- fplot(so_int, reduction="umap", features="nFeature_RNA") + scale_colour_viridis_c(option="inferno") + ggtitle("Feature")
dplot_comb <- dplot_1 + dplot_2 + dplot_3 + dplot_4 + dplot_5 + dplot_6 + fplot_1 + fplot_2 + fplot_3 + fplot_4 + fplot_5 + plot_layout(ncol=4)
dplot_comb
Scale for 'colour' is already present. Adding another scale for 'colour', which will replace the existing scale. Scale for 'colour' is already present. Adding another scale for 'colour', which will replace the existing scale. Scale for 'colour' is already present. Adding another scale for 'colour', which will replace the existing scale. Scale for 'colour' is already present. Adding another scale for 'colour', which will replace the existing scale. Scale for 'colour' is already present. Adding another scale for 'colour', which will replace the existing scale.
features_8000 <- SelectIntegrationFeatures(object.list=so, nfeatures=8000)
so <- PrepSCTIntegration(object.list=so, anchor.features=features_8000)
anchors <- FindIntegrationAnchors(object.list=so, normalization.method="SCT", anchor.features=features_8000, verbose=FALSE)
so_int <- IntegrateData(anchorset=anchors, normalization.method="SCT", verbose=FALSE)
so_int <- RunPCA(so_int, npcs=100, assay="integrated", verbose=FALSE)
so_int <- FindNeighbors(so_int, dims=1:30, k.param=50, verbose=FALSE)
so_int <- FindClusters(so_int, verbose=FALSE, resolution=1, algorithm=1, group.singletons=TRUE)
so_int <- RunUMAP(so_int, dims=1:100, n.neighbors=50, verbose=FALSE)
saveRDS(so_int, "data/object/so_sct_int_hvg8000.rds")
# so_int <- readRDS("data/object/so_sct_int_hvg8000.rds")
options(repr.plot.width=20, repr.plot.height=15)
dplot_1 <- dplot(so_int, reduction="umap", group_by="seurat_clusters") + ggtitle("Louvain")
dplot_2 <- dplot(so_int, reduction="umap", group_by="treatment") + scale_color_manual(values=color$treatment) + ggtitle("Treatment")
dplot_3 <- dplot(so_int, reduction="umap", group_by="sample_rep") + ggtitle("Replicate")
dplot_4 <- dplot(so_int, reduction="umap", group_by="tissue") + scale_color_manual(values=color$tissue) + ggtitle("Tissue")
dplot_5 <- dplot(so_int, reduction="umap", group_by="cc_phase_class") + scale_color_manual(values=color$cc_phase_class) + ggtitle("Cell cycle")
dplot_6 <- dplot(so_int, reduction="umap", group_by="label_fine_haemosphere") + scale_color_manual(values=color$label_fine_haemosphere) + ggtitle("Haemosphere label")
fplot_1 <- fplot(so_int, reduction="umap", features="pMt_RNA") + scale_colour_viridis_c(option="inferno") + ggtitle("Percentage Mt")
fplot_2 <- fplot(so_int, reduction="umap", features="pHb_RNA") + scale_colour_viridis_c(option="inferno") + ggtitle("Percentage Hb")
fplot_3 <- fplot(so_int, reduction="umap", features="pRb_RNA") + scale_colour_viridis_c(option="inferno") + ggtitle("Percentage Rb")
fplot_4 <- fplot(so_int, reduction="umap", features="nCount_RNA") + scale_colour_viridis_c(option="inferno") + ggtitle("UMI")
fplot_5 <- fplot(so_int, reduction="umap", features="nFeature_RNA") + scale_colour_viridis_c(option="inferno") + ggtitle("Feature")
dplot_comb <- dplot_1 + dplot_2 + dplot_3 + dplot_4 + dplot_5 + dplot_6 + fplot_1 + fplot_2 + fplot_3 + fplot_4 + fplot_5 + plot_layout(ncol=4)
dplot_comb
Scale for 'colour' is already present. Adding another scale for 'colour', which will replace the existing scale. Scale for 'colour' is already present. Adding another scale for 'colour', which will replace the existing scale. Scale for 'colour' is already present. Adding another scale for 'colour', which will replace the existing scale. Scale for 'colour' is already present. Adding another scale for 'colour', which will replace the existing scale. Scale for 'colour' is already present. Adding another scale for 'colour', which will replace the existing scale.
# Store data as h5ad
adata <- import("anndata", as="adata", convert=FALSE)
pd <- import("pandas", as="pd", convert=FALSE)
np <- import("numpy", as="np", convert=FALSE)
# Re-load SCT object so that it is not PrepSCTIntegration modified
so <- readRDS("data/object/so_sct.rds")
# Get SCT scale.data
X_1 <- GetAssayData(so[[1]], assay="SCT", slot="scale.data") %>% as.matrix() %>% t()
X_2 <- GetAssayData(so[[2]], assay="SCT", slot="scale.data") %>% as.matrix() %>% t()
X_3 <- GetAssayData(so[[3]], assay="SCT", slot="scale.data") %>% as.matrix() %>% t()
X_4 <- GetAssayData(so[[4]], assay="SCT", slot="scale.data") %>% as.matrix() %>% t()
features <- c(colnames(X_1), colnames(X_2), colnames(X_3), colnames(X_4))
features <- table(features)
features <- features[features==4] %>% names()
X_1 <- X_1[, features]
X_2 <- X_2[, features]
X_3 <- X_3[, features]
X_4 <- X_4[, features]
X <- rbind(X_1, X_2, X_3, X_4)
obs <- rbind(
so[[1]]@meta.data,
so[[2]]@meta.data,
so[[3]]@meta.data,
so[[4]]@meta.data
)
uns <- list(
hvg_int_2000 = features_2000,
hvg_int_4000 = features_4000,
hvg_int_6000 = features_6000,
hvg_int_8000 = features_8000
)
adata_sct <- adata$AnnData(X=X, obs=obs)
adata_sct$var_names <- colnames(X)
adata_sct$uns <- uns
# Get raw counts
so <- readRDS("data/object/qc.rds")
X_raw <- GetAssayData(so, assay="RNA", slot="counts") %>% as.matrix() %>% t()
adata_raw <- adata$AnnData(X=X_raw, obs=obs)
adata_raw$var_names <- colnames(X_raw)
adata_sct$raw <- adata_raw
adata_sct$write_h5ad("data/object/so_sct.h5ad")
None
so <- readRDS("data/object/qc.rds")
so$integrate <- paste0(so$treatment, "-", so$sample_rep)
so <- SplitObject(so, split.by="integrate")
so <- so[c("NaCl-Rep2", "NaCl-Rep1", "CpG-Rep2", "CpG-Rep1")]
so_reg <- lapply(so, function(so) {
so <- SCTransform(so, return.only.var.genes=FALSE, variable.features.n=8000, vars.to.regress=c("msCC_diff_RNA"), verbose=FALSE)
so <- RunPCA(so, npcs=100, assay="SCT", verbose=FALSE)
so <- FindNeighbors(so, dims=1:30, k.param=50, verbose=FALSE)
so <- FindClusters(so, verbose=FALSE, resolution=1, algorithm=1, group.singletons=TRUE)
so <- RunUMAP(so, dims=1:100, n.neighbors=50, verbose=FALSE)
return(so)
}
)
saveRDS(so_reg, "data/object/so_sct_reg.rds")
# so_reg <- readRDS("data/object/so_sct_reg.rds")
options(repr.plot.width=20, repr.plot.height=15)
dplot_1 <- dplot(so_reg[[1]], reduction="umap", group_by="seurat_clusters") + ggtitle("Louvain")
dplot_2 <- dplot(so_reg[[1]], reduction="umap", group_by="treatment") + scale_color_manual(values=color$treatment) + ggtitle("Treatment")
dplot_3 <- dplot(so_reg[[1]], reduction="umap", group_by="sample_rep") + ggtitle("Replicate")
dplot_4 <- dplot(so_reg[[1]], reduction="umap", group_by="tissue") + scale_color_manual(values=color$tissue) + ggtitle("Tissue")
dplot_5 <- dplot(so_reg[[1]], reduction="umap", group_by="cc_phase_class") + scale_color_manual(values=color$cc_phase_class) + ggtitle("Cell cycle")
dplot_6 <- dplot(so_reg[[1]], reduction="umap", group_by="label_fine_haemosphere") + scale_color_manual(values=color$label_fine_haemosphere) + ggtitle("Haemosphere label")
fplot_1 <- fplot(so_reg[[1]], reduction="umap", features="pMt_RNA") + scale_colour_viridis_c(option="inferno") + ggtitle("Percentage Mt")
fplot_2 <- fplot(so_reg[[1]], reduction="umap", features="pHb_RNA") + scale_colour_viridis_c(option="inferno") + ggtitle("Percentage Hb")
fplot_3 <- fplot(so_reg[[1]], reduction="umap", features="pRb_RNA") + scale_colour_viridis_c(option="inferno") + ggtitle("Percentage Rb")
fplot_4 <- fplot(so_reg[[1]], reduction="umap", features="nCount_RNA") + scale_colour_viridis_c(option="inferno") + ggtitle("UMI")
fplot_5 <- fplot(so_reg[[1]], reduction="umap", features="nFeature_RNA") + scale_colour_viridis_c(option="inferno") + ggtitle("Feature")
dplot_comb <- dplot_1 + dplot_2 + dplot_3 + dplot_4 + dplot_5 + dplot_6 + fplot_1 + fplot_2 + fplot_3 + fplot_4 + fplot_5 + plot_layout(ncol=4)
dplot_comb
Scale for 'colour' is already present. Adding another scale for 'colour', which will replace the existing scale. Scale for 'colour' is already present. Adding another scale for 'colour', which will replace the existing scale. Scale for 'colour' is already present. Adding another scale for 'colour', which will replace the existing scale. Scale for 'colour' is already present. Adding another scale for 'colour', which will replace the existing scale. Scale for 'colour' is already present. Adding another scale for 'colour', which will replace the existing scale.
options(repr.plot.width=20, repr.plot.height=15)
dplot_1 <- dplot(so_reg[[2]], reduction="umap", group_by="seurat_clusters") + ggtitle("Louvain")
dplot_2 <- dplot(so_reg[[2]], reduction="umap", group_by="treatment") + scale_color_manual(values=color$treatment) + ggtitle("Treatment")
dplot_3 <- dplot(so_reg[[2]], reduction="umap", group_by="sample_rep") + ggtitle("Replicate")
dplot_4 <- dplot(so_reg[[2]], reduction="umap", group_by="tissue") + scale_color_manual(values=color$tissue) + ggtitle("Tissue")
dplot_5 <- dplot(so_reg[[2]], reduction="umap", group_by="cc_phase_class") + scale_color_manual(values=color$cc_phase_class) + ggtitle("Cell cycle")
dplot_6 <- dplot(so_reg[[2]], reduction="umap", group_by="label_fine_haemosphere") + scale_color_manual(values=color$label_fine_haemosphere) + ggtitle("Haemosphere label")
fplot_1 <- fplot(so_reg[[2]], reduction="umap", features="pMt_RNA") + scale_colour_viridis_c(option="inferno") + ggtitle("Percentage Mt")
fplot_2 <- fplot(so_reg[[2]], reduction="umap", features="pHb_RNA") + scale_colour_viridis_c(option="inferno") + ggtitle("Percentage Hb")
fplot_3 <- fplot(so_reg[[2]], reduction="umap", features="pRb_RNA") + scale_colour_viridis_c(option="inferno") + ggtitle("Percentage Rb")
fplot_4 <- fplot(so_reg[[2]], reduction="umap", features="nCount_RNA") + scale_colour_viridis_c(option="inferno") + ggtitle("UMI")
fplot_5 <- fplot(so_reg[[2]], reduction="umap", features="nFeature_RNA") + scale_colour_viridis_c(option="inferno") + ggtitle("Feature")
dplot_comb <- dplot_1 + dplot_2 + dplot_3 + dplot_4 + dplot_5 + dplot_6 + fplot_1 + fplot_2 + fplot_3 + fplot_4 + fplot_5 + plot_layout(ncol=4)
dplot_comb
Scale for 'colour' is already present. Adding another scale for 'colour', which will replace the existing scale. Scale for 'colour' is already present. Adding another scale for 'colour', which will replace the existing scale. Scale for 'colour' is already present. Adding another scale for 'colour', which will replace the existing scale. Scale for 'colour' is already present. Adding another scale for 'colour', which will replace the existing scale. Scale for 'colour' is already present. Adding another scale for 'colour', which will replace the existing scale.
options(repr.plot.width=20, repr.plot.height=15)
dplot_1 <- dplot(so_reg[[3]], reduction="umap", group_by="seurat_clusters") + ggtitle("Louvain")
dplot_2 <- dplot(so_reg[[3]], reduction="umap", group_by="treatment") + scale_color_manual(values=color$treatment) + ggtitle("Treatment")
dplot_3 <- dplot(so_reg[[3]], reduction="umap", group_by="sample_rep") + ggtitle("Replicate")
dplot_4 <- dplot(so_reg[[3]], reduction="umap", group_by="tissue") + scale_color_manual(values=color$tissue) + ggtitle("Tissue")
dplot_5 <- dplot(so_reg[[3]], reduction="umap", group_by="cc_phase_class") + scale_color_manual(values=color$cc_phase_class) + ggtitle("Cell cycle")
dplot_6 <- dplot(so_reg[[3]], reduction="umap", group_by="label_fine_haemosphere") + scale_color_manual(values=color$label_fine_haemosphere) + ggtitle("Haemosphere label")
fplot_1 <- fplot(so_reg[[3]], reduction="umap", features="pMt_RNA") + scale_colour_viridis_c(option="inferno") + ggtitle("Percentage Mt")
fplot_2 <- fplot(so_reg[[3]], reduction="umap", features="pHb_RNA") + scale_colour_viridis_c(option="inferno") + ggtitle("Percentage Hb")
fplot_3 <- fplot(so_reg[[3]], reduction="umap", features="pRb_RNA") + scale_colour_viridis_c(option="inferno") + ggtitle("Percentage Rb")
fplot_4 <- fplot(so_reg[[3]], reduction="umap", features="nCount_RNA") + scale_colour_viridis_c(option="inferno") + ggtitle("UMI")
fplot_5 <- fplot(so_reg[[4]], reduction="umap", features="nFeature_RNA") + scale_colour_viridis_c(option="inferno") + ggtitle("Feature")
dplot_comb <- dplot_1 + dplot_2 + dplot_3 + dplot_4 + dplot_5 + dplot_6 + fplot_1 + fplot_2 + fplot_3 + fplot_4 + fplot_5 + plot_layout(ncol=4)
dplot_comb
Scale for 'colour' is already present. Adding another scale for 'colour', which will replace the existing scale. Scale for 'colour' is already present. Adding another scale for 'colour', which will replace the existing scale. Scale for 'colour' is already present. Adding another scale for 'colour', which will replace the existing scale. Scale for 'colour' is already present. Adding another scale for 'colour', which will replace the existing scale. Scale for 'colour' is already present. Adding another scale for 'colour', which will replace the existing scale.
options(repr.plot.width=20, repr.plot.height=15)
dplot_1 <- dplot(so_reg[[4]], reduction="umap", group_by="seurat_clusters") + ggtitle("Louvain")
dplot_2 <- dplot(so_reg[[4]], reduction="umap", group_by="treatment") + scale_color_manual(values=color$treatment) + ggtitle("Treatment")
dplot_3 <- dplot(so_reg[[4]], reduction="umap", group_by="sample_rep") + ggtitle("Replicate")
dplot_4 <- dplot(so_reg[[4]], reduction="umap", group_by="tissue") + scale_color_manual(values=color$tissue) + ggtitle("Tissue")
dplot_5 <- dplot(so_reg[[4]], reduction="umap", group_by="cc_phase_class") + scale_color_manual(values=color$cc_phase_class) + ggtitle("Cell cycle")
dplot_6 <- dplot(so_reg[[4]], reduction="umap", group_by="label_fine_haemosphere") + scale_color_manual(values=color$label_fine_haemosphere) + ggtitle("Haemosphere label")
fplot_1 <- fplot(so_reg[[4]], reduction="umap", features="pMt_RNA") + scale_colour_viridis_c(option="inferno") + ggtitle("Percentage Mt")
fplot_2 <- fplot(so_reg[[4]], reduction="umap", features="pHb_RNA") + scale_colour_viridis_c(option="inferno") + ggtitle("Percentage Hb")
fplot_3 <- fplot(so_reg[[4]], reduction="umap", features="pRb_RNA") + scale_colour_viridis_c(option="inferno") + ggtitle("Percentage Rb")
fplot_4 <- fplot(so_reg[[4]], reduction="umap", features="nCount_RNA") + scale_colour_viridis_c(option="inferno") + ggtitle("UMI")
fplot_5 <- fplot(so_reg[[4]], reduction="umap", features="nFeature_RNA") + scale_colour_viridis_c(option="inferno") + ggtitle("Feature")
dplot_comb <- dplot_1 + dplot_2 + dplot_3 + dplot_4 + dplot_5 + dplot_6 + fplot_1 + fplot_2 + fplot_3 + fplot_4 + fplot_5 + plot_layout(ncol=4)
dplot_comb
Scale for 'colour' is already present. Adding another scale for 'colour', which will replace the existing scale. Scale for 'colour' is already present. Adding another scale for 'colour', which will replace the existing scale. Scale for 'colour' is already present. Adding another scale for 'colour', which will replace the existing scale. Scale for 'colour' is already present. Adding another scale for 'colour', which will replace the existing scale. Scale for 'colour' is already present. Adding another scale for 'colour', which will replace the existing scale.
features_2000_reg <- SelectIntegrationFeatures(object.list=so_reg, nfeatures=2000)
so_reg <- PrepSCTIntegration(object.list=so_reg, anchor.features=features_2000_reg)
anchors_reg <- FindIntegrationAnchors(object.list=so_reg, normalization.method="SCT", anchor.features=features_2000_reg, verbose=FALSE)
so_reg_int <- IntegrateData(anchorset=anchors_reg, normalization.method="SCT", verbose=FALSE)
so_reg_int <- RunPCA(so_reg_int, npcs=100, assay="integrated", verbose=FALSE)
so_reg_int <- FindNeighbors(so_reg_int, dims=1:30, k.param=50, verbose=FALSE)
so_reg_int <- FindClusters(so_reg_int, verbose=FALSE, resolution=1, algorithm=1, group.singletons=TRUE)
so_reg_int <- RunUMAP(so_reg_int, dims=1:100, n.neighbors=50, verbose=FALSE)
saveRDS(so_reg_int, "data/object/so_sct_reg_int_hvg_2000.rds")
# so_reg_int <- readRDS("data/object/so_sct_reg_int_hvg_2000.rds")
options(repr.plot.width=20, repr.plot.height=15)
dplot_1 <- dplot(so_reg_int, reduction="umap", group_by="seurat_clusters") + ggtitle("Louvain")
dplot_2 <- dplot(so_reg_int, reduction="umap", group_by="treatment") + scale_color_manual(values=color$treatment) + ggtitle("Treatment")
dplot_3 <- dplot(so_reg_int, reduction="umap", group_by="sample_rep") + ggtitle("Replicate")
dplot_4 <- dplot(so_reg_int, reduction="umap", group_by="tissue") + scale_color_manual(values=color$tissue) + ggtitle("Tissue")
dplot_5 <- dplot(so_reg_int, reduction="umap", group_by="cc_phase_class") + scale_color_manual(values=color$cc_phase_class) + ggtitle("Cell cycle")
dplot_6 <- dplot(so_reg_int, reduction="umap", group_by="label_fine_haemosphere") + scale_color_manual(values=color$label_fine_haemosphere) + ggtitle("Haemosphere label")
fplot_1 <- fplot(so_reg_int, reduction="umap", features="pMt_RNA") + scale_colour_viridis_c(option="inferno") + ggtitle("Percentage Mt")
fplot_2 <- fplot(so_reg_int, reduction="umap", features="pHb_RNA") + scale_colour_viridis_c(option="inferno") + ggtitle("Percentage Hb")
fplot_3 <- fplot(so_reg_int, reduction="umap", features="pRb_RNA") + scale_colour_viridis_c(option="inferno") + ggtitle("Percentage Rb")
fplot_4 <- fplot(so_reg_int, reduction="umap", features="nCount_RNA") + scale_colour_viridis_c(option="inferno") + ggtitle("UMI")
fplot_5 <- fplot(so_reg_int, reduction="umap", features="nFeature_RNA") + scale_colour_viridis_c(option="inferno") + ggtitle("Feature")
dplot_comb <- dplot_1 + dplot_2 + dplot_3 + dplot_4 + dplot_5 + dplot_6 + fplot_1 + fplot_2 + fplot_3 + fplot_4 + fplot_5 + plot_layout(ncol=4)
dplot_comb
Scale for 'colour' is already present. Adding another scale for 'colour', which will replace the existing scale. Scale for 'colour' is already present. Adding another scale for 'colour', which will replace the existing scale. Scale for 'colour' is already present. Adding another scale for 'colour', which will replace the existing scale. Scale for 'colour' is already present. Adding another scale for 'colour', which will replace the existing scale. Scale for 'colour' is already present. Adding another scale for 'colour', which will replace the existing scale.
features_4000_reg <- SelectIntegrationFeatures(object.list=so_reg, nfeatures=4000)
so_reg <- PrepSCTIntegration(object.list=so_reg, anchor.features=features_4000_reg)
anchors_reg <- FindIntegrationAnchors(object.list=so_reg, normalization.method="SCT", anchor.features=features_4000_reg, verbose=FALSE)
so_reg_int <- IntegrateData(anchorset=anchors_reg, normalization.method="SCT", verbose=FALSE)
so_reg_int <- RunPCA(so_reg_int, npcs=100, assay="integrated", verbose=FALSE)
so_reg_int <- FindNeighbors(so_reg_int, dims=1:30, k.param=50, verbose=FALSE)
so_reg_int <- FindClusters(so_reg_int, verbose=FALSE, resolution=1, algorithm=1, group.singletons=TRUE)
so_reg_int <- RunUMAP(so_reg_int, dims=1:100, n.neighbors=50, verbose=FALSE)
saveRDS(so_reg_int, "data/object/so_sct_reg_int_hvg4000.rds")
# so_reg_int <- readRDS("data/object/so_sct_reg_int_hvg4000.rds")
options(repr.plot.width=20, repr.plot.height=15)
dplot_1 <- dplot(so_reg_int, reduction="umap", group_by="seurat_clusters") + ggtitle("Louvain")
dplot_2 <- dplot(so_reg_int, reduction="umap", group_by="treatment") + scale_color_manual(values=color$treatment) + ggtitle("Treatment")
dplot_3 <- dplot(so_reg_int, reduction="umap", group_by="sample_rep") + ggtitle("Replicate")
dplot_4 <- dplot(so_reg_int, reduction="umap", group_by="tissue") + scale_color_manual(values=color$tissue) + ggtitle("Tissue")
dplot_5 <- dplot(so_reg_int, reduction="umap", group_by="cc_phase_class") + scale_color_manual(values=color$cc_phase_class) + ggtitle("Cell cycle")
dplot_6 <- dplot(so_reg_int, reduction="umap", group_by="label_fine_haemosphere") + scale_color_manual(values=color$label_fine_haemosphere) + ggtitle("Haemosphere label")
fplot_1 <- fplot(so_reg_int, reduction="umap", features="pMt_RNA") + scale_colour_viridis_c(option="inferno") + ggtitle("Percentage Mt")
fplot_2 <- fplot(so_reg_int, reduction="umap", features="pHb_RNA") + scale_colour_viridis_c(option="inferno") + ggtitle("Percentage Hb")
fplot_3 <- fplot(so_reg_int, reduction="umap", features="pRb_RNA") + scale_colour_viridis_c(option="inferno") + ggtitle("Percentage Rb")
fplot_4 <- fplot(so_reg_int, reduction="umap", features="nCount_RNA") + scale_colour_viridis_c(option="inferno") + ggtitle("UMI")
fplot_5 <- fplot(so_reg_int, reduction="umap", features="nFeature_RNA") + scale_colour_viridis_c(option="inferno") + ggtitle("Feature")
dplot_comb <- dplot_1 + dplot_2 + dplot_3 + dplot_4 + dplot_5 + dplot_6 + fplot_1 + fplot_2 + fplot_3 + fplot_4 + fplot_5 + plot_layout(ncol=4)
dplot_comb
Scale for 'colour' is already present. Adding another scale for 'colour', which will replace the existing scale. Scale for 'colour' is already present. Adding another scale for 'colour', which will replace the existing scale. Scale for 'colour' is already present. Adding another scale for 'colour', which will replace the existing scale. Scale for 'colour' is already present. Adding another scale for 'colour', which will replace the existing scale. Scale for 'colour' is already present. Adding another scale for 'colour', which will replace the existing scale.
features_6000_reg <- SelectIntegrationFeatures(object.list=so_reg, nfeatures=6000)
so_reg <- PrepSCTIntegration(object.list=so_reg, anchor.features=features_6000_reg)
anchors_reg <- FindIntegrationAnchors(object.list=so_reg, normalization.method="SCT", anchor.features=features_6000_reg, verbose=FALSE)
so_reg_int <- IntegrateData(anchorset=anchors_reg, normalization.method="SCT", verbose=FALSE)
so_reg_int <- RunPCA(so_reg_int, npcs=100, assay="integrated", verbose=FALSE)
so_reg_int <- FindNeighbors(so_reg_int, dims=1:30, k.param=50, verbose=FALSE)
so_reg_int <- FindClusters(so_reg_int, verbose=FALSE, resolution=1, algorithm=1, group.singletons=TRUE)
so_reg_int <- RunUMAP(so_reg_int, dims=1:100, n.neighbors=50, verbose=FALSE)
saveRDS(so_reg_int, "data/object/so_sct_reg_int_hvg6000.rds")
# so_reg_int <- readRDS("data/object/so_sct_reg_int_hvg6000.rds")
options(repr.plot.width=20, repr.plot.height=15)
dplot_1 <- dplot(so_reg_int, reduction="umap", group_by="seurat_clusters") + ggtitle("Louvain")
dplot_2 <- dplot(so_reg_int, reduction="umap", group_by="treatment") + scale_color_manual(values=color$treatment) + ggtitle("Treatment")
dplot_3 <- dplot(so_reg_int, reduction="umap", group_by="sample_rep") + ggtitle("Replicate")
dplot_4 <- dplot(so_reg_int, reduction="umap", group_by="tissue") + scale_color_manual(values=color$tissue) + ggtitle("Tissue")
dplot_5 <- dplot(so_reg_int, reduction="umap", group_by="cc_phase_class") + scale_color_manual(values=color$cc_phase_class) + ggtitle("Cell cycle")
dplot_6 <- dplot(so_reg_int, reduction="umap", group_by="label_fine_haemosphere") + scale_color_manual(values=color$label_fine_haemosphere) + ggtitle("Haemosphere label")
fplot_1 <- fplot(so_reg_int, reduction="umap", features="pMt_RNA") + scale_colour_viridis_c(option="inferno") + ggtitle("Percentage Mt")
fplot_2 <- fplot(so_reg_int, reduction="umap", features="pHb_RNA") + scale_colour_viridis_c(option="inferno") + ggtitle("Percentage Hb")
fplot_3 <- fplot(so_reg_int, reduction="umap", features="pRb_RNA") + scale_colour_viridis_c(option="inferno") + ggtitle("Percentage Rb")
fplot_4 <- fplot(so_reg_int, reduction="umap", features="nCount_RNA") + scale_colour_viridis_c(option="inferno") + ggtitle("UMI")
fplot_5 <- fplot(so_reg_int, reduction="umap", features="nFeature_RNA") + scale_colour_viridis_c(option="inferno") + ggtitle("Feature")
dplot_comb <- dplot_1 + dplot_2 + dplot_3 + dplot_4 + dplot_5 + dplot_6 + fplot_1 + fplot_2 + fplot_3 + fplot_4 + fplot_5 + plot_layout(ncol=4)
dplot_comb
Scale for 'colour' is already present. Adding another scale for 'colour', which will replace the existing scale. Scale for 'colour' is already present. Adding another scale for 'colour', which will replace the existing scale. Scale for 'colour' is already present. Adding another scale for 'colour', which will replace the existing scale. Scale for 'colour' is already present. Adding another scale for 'colour', which will replace the existing scale. Scale for 'colour' is already present. Adding another scale for 'colour', which will replace the existing scale.
features_8000_reg <- SelectIntegrationFeatures(object.list=so_reg, nfeatures=8000)
so_reg <- PrepSCTIntegration(object.list=so_reg, anchor.features=features_8000_reg)
anchors_reg <- FindIntegrationAnchors(object.list=so_reg, normalization.method="SCT", anchor.features=features_8000_reg, verbose=FALSE)
so_reg_int <- IntegrateData(anchorset=anchors_reg, normalization.method="SCT", verbose=FALSE)
so_reg_int <- RunPCA(so_reg_int, npcs=100, assay="integrated", verbose=FALSE)
so_reg_int <- FindNeighbors(so_reg_int, dims=1:30, k.param=50, verbose=FALSE)
so_reg_int <- FindClusters(so_reg_int, verbose=FALSE, resolution=1, algorithm=1, group.singletons=TRUE)
so_reg_int <- RunUMAP(so_reg_int, dims=1:100, n.neighbors=50, verbose=FALSE)
saveRDS(so_reg_int, "data/object/so_sct_reg_int_hvg8000.rds")
# so_reg_int <- readRDS("data/object/so_sct_reg_int_hvg8000.rds")
options(repr.plot.width=20, repr.plot.height=15)
dplot_1 <- dplot(so_reg_int, reduction="umap", group_by="seurat_clusters") + ggtitle("Louvain")
dplot_2 <- dplot(so_reg_int, reduction="umap", group_by="treatment") + scale_color_manual(values=color$treatment) + ggtitle("Treatment")
dplot_3 <- dplot(so_reg_int, reduction="umap", group_by="sample_rep") + ggtitle("Replicate")
dplot_4 <- dplot(so_reg_int, reduction="umap", group_by="tissue") + scale_color_manual(values=color$tissue) + ggtitle("Tissue")
dplot_5 <- dplot(so_reg_int, reduction="umap", group_by="cc_phase_class") + scale_color_manual(values=color$cc_phase_class) + ggtitle("Cell cycle")
dplot_6 <- dplot(so_reg_int, reduction="umap", group_by="label_fine_haemosphere") + scale_color_manual(values=color$label_fine_haemosphere) + ggtitle("Haemosphere label")
fplot_1 <- fplot(so_reg_int, reduction="umap", features="pMt_RNA") + scale_colour_viridis_c(option="inferno") + ggtitle("Percentage Mt")
fplot_2 <- fplot(so_reg_int, reduction="umap", features="pHb_RNA") + scale_colour_viridis_c(option="inferno") + ggtitle("Percentage Hb")
fplot_3 <- fplot(so_reg_int, reduction="umap", features="pRb_RNA") + scale_colour_viridis_c(option="inferno") + ggtitle("Percentage Rb")
fplot_4 <- fplot(so_reg_int, reduction="umap", features="nCount_RNA") + scale_colour_viridis_c(option="inferno") + ggtitle("UMI")
fplot_5 <- fplot(so_reg_int, reduction="umap", features="nFeature_RNA") + scale_colour_viridis_c(option="inferno") + ggtitle("Feature")
dplot_comb <- dplot_1 + dplot_2 + dplot_3 + dplot_4 + dplot_5 + dplot_6 + fplot_1 + fplot_2 + fplot_3 + fplot_4 + fplot_5 + plot_layout(ncol=4)
dplot_comb
Scale for 'colour' is already present. Adding another scale for 'colour', which will replace the existing scale. Scale for 'colour' is already present. Adding another scale for 'colour', which will replace the existing scale. Scale for 'colour' is already present. Adding another scale for 'colour', which will replace the existing scale. Scale for 'colour' is already present. Adding another scale for 'colour', which will replace the existing scale. Scale for 'colour' is already present. Adding another scale for 'colour', which will replace the existing scale.
# Store data as h5ad
adata <- import("anndata", as="adata", convert=FALSE)
pd <- import("pandas", as="pd", convert=FALSE)
np <- import("numpy", as="np", convert=FALSE)
# Re-load SCT object so that it is not PrepSCTIntegration modified
so <- readRDS("data/object/so_sct_reg.rds")
# Get SCT scale.data
X_1 <- GetAssayData(so[[1]], assay="SCT", slot="scale.data") %>% as.matrix() %>% t()
X_2 <- GetAssayData(so[[2]], assay="SCT", slot="scale.data") %>% as.matrix() %>% t()
X_3 <- GetAssayData(so[[3]], assay="SCT", slot="scale.data") %>% as.matrix() %>% t()
X_4 <- GetAssayData(so[[4]], assay="SCT", slot="scale.data") %>% as.matrix() %>% t()
features <- c(colnames(X_1), colnames(X_2), colnames(X_3), colnames(X_4))
features <- table(features)
features <- features[features==4] %>% names()
X_1 <- X_1[, features]
X_2 <- X_2[, features]
X_3 <- X_3[, features]
X_4 <- X_4[, features]
X <- rbind(X_1, X_2, X_3, X_4)
obs <- rbind(
so[[1]]@meta.data,
so[[2]]@meta.data,
so[[3]]@meta.data,
so[[4]]@meta.data
)
uns <- list(
hvg_int_2000 = features_2000_reg,
hvg_int_4000 = features_4000_reg,
hvg_int_6000 = features_6000_reg,
hvg_int_8000 = features_8000_reg
)
adata_sct <- adata$AnnData(X=X, obs=obs)
adata_sct$var_names <- colnames(X)
adata_sct$uns <- uns
# Get raw counts
so <- readRDS("data/object/qc.rds")
X_raw <- GetAssayData(so, assay="RNA", slot="counts") %>% as.matrix() %>% t()
adata_raw <- adata$AnnData(X=X_raw, obs=obs)
adata_raw$var_names <- colnames(X_raw)
adata_sct$raw <- adata_raw
adata_sct$write_h5ad("data/object/so_sct_reg.h5ad")
None